home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1993 July / InfoMagic USENET CD-ROM July 1993.ISO / sources / unix / volume18 / mtvraytrace / part02 < prev    next >
Encoding:
Internet Message Format  |  1989-03-26  |  52.7 KB
  1. Subject:  v18i071:  A ray-tracing package, Part02/03
  2. Newsgroups: comp.sources.unix
  3. Sender: sources
  4. Approved: rsalz@uunet.UU.NET
  5.  
  6. Submitted-by: Mark VandeWettering <markv@drizzle.cs.uoregon.edu>
  7. Posting-number: Volume 18, Issue 71
  8. Archive-name: mtvraytrace/part02
  9.  
  10. #! /bin/sh
  11. # This is a shell archive.  Remove anything before this line, then unpack
  12. # it by saving it into a file and typing "sh file".  To overwrite existing
  13. # files, type "sh file -c".  You can also feed this as standard input via
  14. # unshar, or by typing "sh <file", e.g..  If this archive is complete, you
  15. # will see the following message at the end:
  16. #        "End of archive 2 (of 3)."
  17. # Contents:  NFF color.c cone.c intersect.c poly.c screen.c shade.c
  18. #   tri.c
  19. # Wrapped by markv@tillamook on Mon Nov 14 12:50:26 1988
  20. PATH=/bin:/usr/bin:/usr/ucb ; export PATH
  21. if test -f 'NFF' -a "${1}" != "-c" ; then 
  22.   echo shar: Will not clobber existing file \"'NFF'\"
  23. else
  24. echo shar: Extracting \"'NFF'\" \(8711 characters\)
  25. sed "s/^X//" >'NFF' <<'END_OF_FILE'
  26. XFrom uoregon!uw-beaver!cornell!batcomputer!saponara Mon Oct  3 13:30:35 PDT 1988
  27. XArticle 3096 of comp.graphics:
  28. XPath: uoregon!uw-beaver!cornell!batcomputer!saponara
  29. X>From: saponara@batcomputer.tn.cornell.edu (John Saponara)
  30. XNewsgroups: comp.graphics
  31. XSubject: A description of NFF (Neutral File Format)
  32. XSummary: here you go...
  33. XKeywords: NFF, benchmark
  34. XMessage-ID: <6443@batcomputer.tn.cornell.edu>
  35. XDate: 3 Oct 88 19:57:58 GMT
  36. XReply-To: saponara@tcgould.tn.cornell.edu (Eric Haines, in truth)
  37. XOrganization: 3D/Eye Inc
  38. XLines: 235
  39. XReferences:
  40. XStatus: OR
  41. X
  42. XHere's my first draft, hot off the disk.
  43. X
  44. X    Eric Haines (not John Saponara)
  45. X
  46. XNote: don't reply to me here, as this account disappears soon.  Instead,
  47. X    write me at hpfcla!hpfcrs!eye!erich@hplabs.hp.com
  48. X
  49. X-----------------------------------------------------------------------------
  50. X
  51. XNeutral File Format (NFF),
  52. X    by Eric Haines, 3D/Eye Inc, 410 E Upland Rd, Ithaca, NY 14850
  53. X    (607)-257-1381
  54. X    email: hpfcla!hpfcrs!eye!erich@hplabs.hp.com
  55. X
  56. X
  57. X[This is a description of the format used in the SPD package.  Any comments
  58. Xon how to expand this format are appreciated.  Some extensions seem obvious
  59. Xto me (e.g. adding directional lights, circles, and tori), but I want to take
  60. Xmy time, gather opinions, and get it more-or-less right the first go round.
  61. XAlso, I'm ridiculously busy right now. -EAH]
  62. X
  63. XDraft document #1, 10/3/88
  64. X
  65. XThe NFF (Neutral File Format) is designed as a minimal scene description
  66. Xlanguage.  The language was designed in order to test various rendering
  67. Xalgorithms and efficiency schemes.  It is meant to describe the geometry and
  68. Xbasic surface characteristics of objects, the placement of lights, and the
  69. Xviewing frustum for the eye.  Some additional information is provided for
  70. Xesthetic reasons (such as the color of the objects, which is not strictly
  71. Xnecessary for testing rendering algorithms).
  72. X
  73. XFuture enhancements include:  circle and torus objects, spline surfaces
  74. Xwith trimming curves, directional lights, characteristics for positional
  75. Xlights, CSG descriptions, and probably more by the time you read this.
  76. XComments, suggestions, and criticisms are all welcome.
  77. X
  78. XAt present the NFF file format is used in conjunction with the SPD (Standard
  79. XProcedural Database) software, a package designed to create a variety of
  80. Xdatabases for testing rendering schemes.  The SPD package is available
  81. Xfrom Netlib and via ftp from drizzle.cs.uoregon.edu.  For more information
  82. Xabout SPD see "A Proposal for Standard Graphics Environments," IEEE Computer
  83. XGraphics and Applications, vol. 7, no. 11, November 1987, pp. 3-5.
  84. X
  85. XBy providing a minimal interface, NFF is meant to act as a simple format to
  86. Xallow the programmer to quickly write filters to move from NFF to the
  87. Xlocal file format.  Presently the following entities are supported:
  88. X     A simple perspective frustum
  89. X     A positional (vs. directional) light source description
  90. X     A background color description
  91. X     A surface properties description
  92. X     Polygon, polygonal patch, cylinder/cone, and sphere descriptions
  93. X
  94. XFiles are output as lines of text.  For each entity, the first line
  95. Xdefines its type.  The rest of the first line and possibly other lines
  96. Xcontain further information about the entity.  Entities include:
  97. X
  98. X"v"  - viewing vectors and angles
  99. X"l"  - positional light location
  100. X"b"  - background color
  101. X"f"  - object material properties
  102. X"c"  - cone or cylinder primitive
  103. X"s"  - sphere primitive
  104. X"p"  - polygon primitive
  105. X"pp" - polygonal patch primitive
  106. X
  107. X
  108. XThese are explained in depth below:
  109. X
  110. XViewpoint location.  Description:
  111. X    "v"
  112. X    "from" Fx Fy Fz
  113. X    "at" Ax Ay Az
  114. X    "up" Ux Uy Uz
  115. X    "angle" angle
  116. X    "hither" hither
  117. X    "resolution" xres yres
  118. X
  119. XFormat:
  120. X
  121. X    v
  122. X    from %g %g %g
  123. X    at %g %g %g
  124. X    up %g %g %g
  125. X    angle %g
  126. X    hither %g
  127. X    resolution %d %d
  128. X
  129. XThe parameters are:
  130. X
  131. X    From:  the eye location in XYZ.
  132. X    At:    a position to be at the center of the image, in XYZ world
  133. X       coordinates.  A.k.a. "lookat".
  134. X    Up:    a vector defining which direction is up, as an XYZ vector.
  135. X    Angle: in degrees, defined as from the center of top pixel row to
  136. X       bottom pixel row and left column to right column.
  137. X    Resolution: in pixels, in x and in y.
  138. X
  139. X  Note that no assumptions are made about normalizing the data (e.g. the
  140. X  from-at distance does not have to be 1).  Also, vectors are not
  141. X  required to be perpendicular to each other.
  142. X
  143. X  For all databases some viewing parameters are always the same:
  144. X    Yon is "at infinity."
  145. X    Aspect ratio is 1.0.
  146. X
  147. X  A view entity must be defined before any objects are defined (this
  148. X  requirement is so that NFF files can be used by hidden surface machines).
  149. X
  150. X--------
  151. X
  152. XPositional light.  A light is defined by XYZ position.  Description:
  153. X    "b" X Y Z
  154. X
  155. XFormat:
  156. X    l %g %g %g
  157. X
  158. X    All light entities must be defined before any objects are defined (this
  159. X    requirement is so that NFF files can be used by hidden surface machines).
  160. X    Lights have a non-zero intensity of no particular value [this definition
  161. X    may change soon, with the addition of an intensity and/or color].
  162. X
  163. X--------
  164. X
  165. XBackground color.  A color is simply RGB with values between 0 and 1:
  166. X    "b" R G B
  167. X
  168. XFormat:
  169. X    b %g %g %g
  170. X
  171. X    If no background color is set, assume RGB = {0,0,0}.
  172. X
  173. X--------
  174. X
  175. XFill color and shading parameters.  Description:
  176. X     "f" red green blue Kd Ks Shine T index_of_refraction
  177. X
  178. XFormat:
  179. X    f %g %g %g %g %g %g %g %g
  180. X
  181. X    RGB is in terms of 0.0 to 1.0.
  182. X
  183. X    Kd is the diffuse component, Ks the specular, Shine is the Phong cosine
  184. X    power for highlights, T is transmittance (fraction of light passed per
  185. X    unit).  Usually, 0 <= Kd <= 1 and 0 <= Ks <= 1, though it is not required
  186. X    that Kd + Ks == 1.  Note that transmitting objects ( T > 0 ) are considered
  187. X    to have two sides for algorithms that need these (normally objects have
  188. X    one side).
  189. X  
  190. X    The fill color is used to color the objects following it until a new color
  191. X    is assigned.
  192. X
  193. X--------
  194. X
  195. XObjects:  all objects are considered one-sided, unless the second side is
  196. Xneeded for transmittance calculations (e.g. you cannot throw out the second
  197. Xintersection of a transparent sphere in ray tracing).
  198. X
  199. XCylinder or cone.  A cylinder is defined as having a radius and an axis
  200. X    defined by two points, which also define the top and bottom edge of the
  201. X    cylinder.  A cone is defined similarly, the difference being that the apex
  202. X    and base radii are different.  The apex radius is defined as being smaller
  203. X    than the base radius.  Note that the surface exists without endcaps.  The
  204. X    cone or cylinder description:
  205. X
  206. X    "c"
  207. X    base.x base.y base.z base_radius
  208. X    apex.x apex.y apex.z apex_radius
  209. X
  210. XFormat:
  211. X    c
  212. X    %g %g %g %g
  213. X    %g %g %g %g
  214. X
  215. X    A negative value for both radii means that only the inside of the object is
  216. X    visible (objects are normally considered one sided, with the outside
  217. X    visible).  Note that the base and apex cannot be coincident for a cylinder
  218. X    or cone.
  219. X
  220. X--------
  221. X
  222. XSphere.  A sphere is defined by a radius and center position:
  223. X    "s" center.x center.y center.z radius
  224. X
  225. XFormat:
  226. X    s %g %g %g %g
  227. X
  228. X    If the radius is negative, then only the sphere's inside is visible
  229. X    (objects are normally considered one sided, with the outside visible).
  230. X
  231. X--------
  232. X
  233. XPolygon.  A polygon is defined by a set of vertices.  With these databases,
  234. X    a polygon is defined to have all points coplanar.  A polygon has only
  235. X    one side, with the order of the vertices being counterclockwise as you
  236. X    face the polygon (right-handed coordinate system).  The first two edges
  237. X    must form a non-zero convex angle, so that the normal and side visibility
  238. X    can be determined.  Description:
  239. X
  240. X    "p" total_vertices
  241. X    vert1.x vert1.y vert1.z
  242. X    [etc. for total_vertices vertices]
  243. X
  244. XFormat:
  245. X    p %d
  246. X    [ %g %g %g ] <-- for total_vertices vertices
  247. X
  248. X--------
  249. X
  250. XPolygonal patch.  A patch is defined by a set of vertices and their normals.
  251. X    With these databases, a patch is defined to have all points coplanar.
  252. X    A patch has only one side, with the order of the vertices being
  253. X    counterclockwise as you face the patch (right-handed coordinate system).
  254. X    The first two edges must form a non-zero convex angle, so that the normal
  255. X    and side visibility can be determined.  Description:
  256. X
  257. X    "pp" total_vertices
  258. X    vert1.x vert1.y vert1.z norm1.x norm1.y norm1.z
  259. X    [etc. for total_vertices vertices]
  260. X
  261. XFormat:
  262. X    pp %d
  263. X    [ %g %g %g %g %g %g ] <-- for total_vertices vertices
  264. X
  265. X--------
  266. X
  267. XComment.  Description:
  268. X    "#" [ string ]
  269. X
  270. XFormat:
  271. X    # [ string ]
  272. X
  273. X    As soon as a "#" character is detected, the rest of the line is considered
  274. X    a comment.
  275. X
  276. X-------------------------------------------------------------------------------
  277. X
  278. X
  279. END_OF_FILE
  280. if test 8711 -ne `wc -c <'NFF'`; then
  281.     echo shar: \"'NFF'\" unpacked with wrong size!
  282. fi
  283. # end of 'NFF'
  284. fi
  285. if test -f 'color.c' -a "${1}" != "-c" ; then 
  286.   echo shar: Will not clobber existing file \"'color.c'\"
  287. else
  288. echo shar: Extracting \"'color.c'\" \(6975 characters\)
  289. sed "s/^X//" >'color.c' <<'END_OF_FILE'
  290. X/***********************************************************************
  291. X * $Author: markv $
  292. X * $Revision: 1.2 $
  293. X * $Date: 88/09/12 12:53:47 $
  294. X * $Log:    color.c,v $
  295. X * Revision 1.2  88/09/12  12:53:47  markv
  296. X * Fixed problem in LookupColorbyName, had return ; and return(0).
  297. X * [ Thank you lint! ]
  298. X * 
  299. X * Revision 1.1  88/09/11  11:00:37  markv
  300. X * Initial revision
  301. X * 
  302. X ***********************************************************************/
  303. X
  304. X#include <stdio.h>
  305. X#include "defs.h"
  306. X
  307. X#define        NCOLORS        (142)
  308. X
  309. Xtypedef struct t_color_entry {
  310. X    char *    ce_name ;
  311. X    Vec     ce_color ;
  312. X} ColorEntry ;
  313. X
  314. X#define LESS_THAN -1
  315. X#define GREATER_THAN 1
  316. X#define EQUAL_TO 0
  317. X
  318. X/*
  319. X * Note: These colors must be in sorted order, because we binary search
  320. X * for them.
  321. X *
  322. X * They were swiped from the X-11 distribution.  Sorry....
  323. X */
  324. X
  325. XColorEntry Colors[] = {
  326. X    "Aquamarine", {.439216, .858824, .576471},
  327. X    "Black", {0, 0, 0},
  328. X    "Blue", {0, 0, 1},
  329. X    "BlueViolet", {.623529, .372549, .623529},
  330. X    "Brown", {.647059, .164706, .164706},
  331. X    "CadetBlue", {.372549, .623529, .623529},
  332. X    "Coral", {1, .498039, 0},
  333. X    "CornflowerBlue", {.258824, .258824, .435294},
  334. X    "Cyan", {0, 1, 1},
  335. X    "DarkGreen", {.184314, .309804, .184314},
  336. X    "DarkOliveGreen", {.309804, .309804, .184314},
  337. X    "DarkOrchid", {.6, .196078, .8},
  338. X    "DarkSlateBlue", {.419608, .137255, .556863},
  339. X    "DarkSlateGray", {.184314, .309804, .309804},
  340. X    "DarkSlateGrey", {.184314, .309804, .309804},
  341. X    "DarkTurquoise", {.439216, .576471, .858824},
  342. X    "DimGray", {.329412, .329412, .329412},
  343. X    "DimGrey", {.329412, .329412, .329412},
  344. X    "Firebrick", {.556863, .137255, .137255},
  345. X    "ForestGreen", {.137255, .556863, .137255},
  346. X    "Gold", {.8, .498039, .196078},
  347. X    "Goldenrod", {.858824, .858824, .439216},
  348. X    "Gray", {.752941, .752941, .752941},
  349. X    "Green", {0, 1, 0},
  350. X    "GreenYellow", {.576471, .858824, .439216},
  351. X    "Grey", {.752941, .752941, .752941},
  352. X    "IndianRed", {.309804, .184314, .184314},
  353. X    "Khaki", {.623529, .623529, .372549},
  354. X    "LightBlue", {.74902, .847059, .847059},
  355. X    "LightGray", {.658824, .658824, .658824},
  356. X    "LightGrey", {.658824, .658824, .658824},
  357. X    "LightSteelBlue", {.560784, .560784, .737255},
  358. X    "LimeGreen", {.196078, .8, .196078},
  359. X    "Magenta", {1, 0, 1},
  360. X    "Maroon", {.556863, .137255, .419608},
  361. X    "MediumAquamarine", {.196078, .8, .6},
  362. X    "MediumBlue", {.196078, .196078, .8},
  363. X    "MediumForestGreen", {.419608, .556863, .137255},
  364. X    "MediumGoldenrod", {.917647, .917647, .678431},
  365. X    "MediumOrchid", {.576471, .439216, .858824},
  366. X    "MediumSeaGreen", {.258824, .435294, .258824},
  367. X    "MediumSlateBlue", {.498039, 0, 1},
  368. X    "MediumSpringGreen", {.498039, 1, 0},
  369. X    "MediumTurquoise", {.439216, .858824, .858824},
  370. X    "MediumVioletRed", {.858824, .439216, .576471},
  371. X    "MidnightBlue", {.184314, .184314, .309804},
  372. X    "Navy", {.137255, .137255, .556863},
  373. X    "NavyBlue", {.137255, .137255, .556863},
  374. X    "Orange", {.8, .196078, .196078},
  375. X    "OrangeRed", {1, 0, .498039},
  376. X    "Orchid", {.858824, .439216, .858824},
  377. X    "PaleGreen", {.560784, .737255, .560784},
  378. X    "Pink", {.737255, .560784, .560784},
  379. X    "Plum", {.917647, .678431, .917647},
  380. X    "Red", {1, 0, 0},
  381. X    "Salmon", {.435294, .258824, .258824},
  382. X    "SeaGreen", {.137255, .556863, .419608},
  383. X    "Sienna", {.556863, .419608, .137255},
  384. X    "SkyBlue", {.196078, .6, .8},
  385. X    "SlateBlue", {0, .498039, 1},
  386. X    "SpringGreen", {0, 1, .498039},
  387. X    "SteelBlue", {.137255, .419608, .556863},
  388. X    "Tan", {.858824, .576471, .439216},
  389. X    "Thistle", {.847059, .74902, .847059},
  390. X    "Turquoise", {.678431, .917647, .917647},
  391. X    "Violet", {.309804, .184314, .309804},
  392. X    "VioletRed", {.8, .196078, .6},
  393. X    "Wheat", {.847059, .847059, .74902},
  394. X    "White", {.988235, .988235, .988235},
  395. X    "Yellow", {1, 1, 0},
  396. X    "YellowGreen", {.6, .8, .196078},
  397. X    "aquamarine", {.439216, .858824, .576471},
  398. X    "black", {0, 0, 0},
  399. X    "blue", {0, 0, 1},
  400. X    "blue_violet", {.623529, .372549, .623529},
  401. X    "brown", {.647059, .164706, .164706},
  402. X    "cadet_blue", {.372549, .623529, .623529},
  403. X    "coral", {1, .498039, 0},
  404. X    "cornflower_blue", {.258824, .258824, .435294},
  405. X    "cyan", {0, 1, 1},
  406. X    "dark_green", {.184314, .309804, .184314},
  407. X    "dark_olive_green", {.309804, .309804, .184314},
  408. X    "dark_orchid", {.6, .196078, .8},
  409. X    "dark_slate_blue", {.419608, .137255, .556863},
  410. X    "dark_slate_gray", {.184314, .309804, .309804},
  411. X    "dark_slate_grey", {.184314, .309804, .309804},
  412. X    "dark_turquoise", {.439216, .576471, .858824},
  413. X    "dim_gray", {.329412, .329412, .329412},
  414. X    "dim_grey", {.329412, .329412, .329412},
  415. X    "firebrick", {.556863, .137255, .137255},
  416. X    "forest_green", {.137255, .556863, .137255},
  417. X    "gold", {.8, .498039, .196078},
  418. X    "goldenrod", {.858824, .858824, .439216},
  419. X    "gray", {.752941, .752941, .752941},
  420. X    "green", {0, 1, 0},
  421. X    "green_yellow", {.576471, .858824, .439216},
  422. X    "grey", {.752941, .752941, .752941},
  423. X    "indian_red", {.309804, .184314, .184314},
  424. X    "khaki", {.623529, .623529, .372549},
  425. X    "light_blue", {.74902, .847059, .847059},
  426. X    "light_gray", {.658824, .658824, .658824},
  427. X    "light_grey", {.658824, .658824, .658824},
  428. X    "light_steel_blue", {.560784, .560784, .737255},
  429. X    "lime_green", {.196078, .8, .196078},
  430. X    "magenta", {1, 0, 1},
  431. X    "maroon", {.556863, .137255, .419608},
  432. X    "medium_aquamarine", {.196078, .8, .6},
  433. X    "medium_blue", {.196078, .196078, .8},
  434. X    "medium_forest_green", {.419608, .556863, .137255},
  435. X    "medium_goldenrod", {.917647, .917647, .678431},
  436. X    "medium_orchid", {.576471, .439216, .858824},
  437. X    "medium_sea_green", {.258824, .435294, .258824},
  438. X    "medium_slate_blue", {.498039, 0, 1},
  439. X    "medium_spring_green", {.498039, 1, 0},
  440. X    "medium_turquoise", {.439216, .858824, .858824},
  441. X    "medium_violet_red", {.858824, .439216, .576471},
  442. X    "midnight_blue", {.184314, .184314, .309804},
  443. X    "navy", {.137255, .137255, .556863},
  444. X    "navy_blue", {.137255, .137255, .556863},
  445. X    "orange", {.8, .196078, .196078},
  446. X    "orange_red", {1, 0, .498039},
  447. X    "orchid", {.858824, .439216, .858824},
  448. X    "pale_green", {.560784, .737255, .560784},
  449. X    "pink", {.737255, .560784, .560784},
  450. X    "plum", {.917647, .678431, .917647},
  451. X    "red", {1, 0, 0},
  452. X    "salmon", {.435294, .258824, .258824},
  453. X    "sea_green", {.137255, .556863, .419608},
  454. X    "sienna", {.556863, .419608, .137255},
  455. X    "sky_blue", {.196078, .6, .8},
  456. X    "slate_blue", {0, .498039, 1},
  457. X    "spring_green", {0, 1, .498039},
  458. X    "steel_blue", {.137255, .419608, .556863},
  459. X    "tan", {.858824, .576471, .439216},
  460. X    "thistle", {.847059, .74902, .847059},
  461. X    "turquoise", {.678431, .917647, .917647},
  462. X    "violet", {.309804, .184314, .309804},
  463. X    "violet_red", {.8, .196078, .6},
  464. X    "wheat", {.847059, .847059, .74902},
  465. X    "white", {.988235, .988235, .988235},
  466. X    "yellow", {1, 1, 0},
  467. X    "yellow_green", {.6, .8, .196078}
  468. X} ;
  469. X
  470. Xint
  471. XLookupColorByName(name, color)
  472. X char * name ;
  473. X Vec color ;
  474. X{
  475. X    int rc ;
  476. X    rc = BinarySearch(name, 0, NCOLORS - 1 , Colors) ;
  477. X    if (rc < 0) {
  478. X        return(0) ;
  479. X    }
  480. X
  481. X    VecCopy(Colors[rc].ce_color, color) ;
  482. X    return 1 ;
  483. X}
  484. X
  485. X
  486. Xint 
  487. XBinarySearch(name, l, h, array)
  488. X char * name ;
  489. X int l, h ;
  490. X ColorEntry array[] ;
  491. X{
  492. X    int m, rc ;
  493. X    if (l > h)
  494. X        return(-1) ;
  495. X    
  496. X    m = (l + h) / 2 ;
  497. X
  498. X    rc = strcmp(name, array[m].ce_name) ;
  499. X    if (rc == 0)
  500. X        return m ;
  501. X    else if (rc < 0)
  502. X        return BinarySearch(name, l, m-1, array) ;
  503. X    else
  504. X        return BinarySearch(name, m + 1, h, array) ;
  505. X}
  506. END_OF_FILE
  507. if test 6975 -ne `wc -c <'color.c'`; then
  508.     echo shar: \"'color.c'\" unpacked with wrong size!
  509. fi
  510. # end of 'color.c'
  511. fi
  512. if test -f 'cone.c' -a "${1}" != "-c" ; then 
  513.   echo shar: Will not clobber existing file \"'cone.c'\"
  514. else
  515. echo shar: Extracting \"'cone.c'\" \(6358 characters\)
  516. sed "s/^X//" >'cone.c' <<'END_OF_FILE'
  517. X/***********************************************************************
  518. X * $Author: markv $
  519. X * $Revision: 1.3 $
  520. X * $Date: 88/10/04 14:30:02 $
  521. X * $Log:    cone.c,v $
  522. X * Revision 1.3  88/10/04  14:30:02  markv
  523. X * Fixed bug reported by koblas@mips.
  524. X * 
  525. X * Revision 1.2  88/09/15  09:33:02  markv
  526. X * Fixed bug reported by koblas@mips.  Cones which were specified with
  527. X * axis coincident with -y axis had problems.  Caused by incorrect 
  528. X * handling when finding orthogonal axes for the local coordinate system.
  529. X * 
  530. X * Revision 1.1  88/09/11  11:00:38  markv
  531. X * Initial revision
  532. X * 
  533. X ***********************************************************************/
  534. X#include <stdio.h>
  535. X#include <math.h>
  536. X#include "defs.h"
  537. X#include "extern.h"
  538. X
  539. Xtypedef struct t_conedata {
  540. X    Vec         cone_base ;
  541. X    Flt        cone_base_radius ;
  542. X    Flt        cone_base_d ;
  543. X    Vec        cone_apex ;
  544. X    Flt        cone_apex_radius ;
  545. X    Vec        cone_u ;
  546. X    Vec        cone_v ;
  547. X    Vec        cone_w ;
  548. X    Flt        cone_height ;
  549. X    Flt        cone_slope ;
  550. X    Flt        cone_min_d ;
  551. X    Flt        cone_max_d ;
  552. X} ConeData ;
  553. X
  554. Xint ConePrint ();
  555. Xint ConeIntersect ();
  556. Xint ConeNormal ();
  557. X
  558. XObjectProcs ConeProcs = {
  559. X    ConePrint,
  560. X    ConeIntersect,
  561. X    ConeNormal,
  562. X} ;
  563. X
  564. Xint 
  565. XConePrint(obj)
  566. X Object *obj ;
  567. X{
  568. X    ConeData * cp ;
  569. X
  570. X    cp = (ConeData *) obj -> o_data ;
  571. X
  572. X    printf("c %g %g %g %g %g %g %g %g\n", cp -> cone_base[0], 
  573. X                  cp -> cone_base[1],
  574. X                  cp -> cone_base[2],
  575. X                  cp -> cone_base_radius,
  576. X                  cp -> cone_apex[0],
  577. X                  cp -> cone_apex[1],
  578. X                  cp -> cone_apex[2],
  579. X                  cp -> cone_apex_radius) ;
  580. X}
  581. X
  582. XConeIntersect(obj, ray, hit)
  583. X Object * obj ;
  584. X Ray * ray ;
  585. X Isect * hit ;
  586. X{
  587. X    Ray tray ;
  588. X    ConeData * cd ;
  589. X    Vec V, P ;
  590. X    Flt a, b, c, d, disc ;
  591. X    Flt t1, t2 ;
  592. X    int nroots ;
  593. X
  594. X    cd = (ConeData *) (obj -> o_data) ;
  595. X
  596. X    /*
  597. X     * First, we get the coordinates of the ray origin in 
  598. X     * the objects space....
  599. X     */
  600. X    
  601. X    VecSub(ray -> P, cd -> cone_base, V) ;
  602. X
  603. X    tray.P[0] = VecDot(V, cd -> cone_u) ;
  604. X    tray.P[1] = VecDot(V, cd -> cone_v) ;
  605. X    tray.P[2] = VecDot(V, cd -> cone_w) ;
  606. X
  607. X    /*
  608. X    VecAdd(ray -> P, ray -> D, V) ;
  609. X    VecSub(V, cd -> cone_base, V) ;
  610. X    */
  611. X
  612. X    tray.D[0] = VecDot(ray -> D, cd -> cone_u) ;
  613. X    tray.D[1] = VecDot(ray -> D, cd -> cone_v) ;
  614. X    tray.D[2] = VecDot(ray -> D, cd -> cone_w) ;
  615. X
  616. X    /*
  617. X    VecSub(tray.D, tray.P, tray.D) ;
  618. X    */
  619. X
  620. X    a = tray.D[0] * tray.D[0] 
  621. X        + tray.D[1] * tray.D[1] 
  622. X        - cd -> cone_slope * cd -> cone_slope * tray.D[2] * tray.D[2] ;
  623. X
  624. X    b = 2.0 * (tray.P[0] * tray.D[0] + tray.P[1] * tray.D[1] - 
  625. X        cd -> cone_slope * cd -> cone_slope * tray.P[2] * tray.D[2]
  626. X        - cd -> cone_base_radius * cd -> cone_slope * tray.D[2]) ;
  627. X
  628. X    c = cd -> cone_slope * tray.P[2] + cd -> cone_base_radius ;
  629. X    c = tray.P[0] * tray.P[0] + tray.P[1] * tray.P[1] - (c * c) ;
  630. X
  631. X    disc = b * b - 4.0 * a * c ;
  632. X
  633. X    if (disc < 0.0)
  634. X        return (0) ;
  635. X    
  636. X    disc = (Flt) sqrt(disc) ;
  637. X    t1 = (-b - disc) / (2.0 * a) ;
  638. X    t2 = (-b + disc) / (2.0 * a) ;
  639. X
  640. X    if (t2 < rayeps)
  641. X        return (0) ;
  642. X    if (t1 < rayeps) {
  643. X        nroots = 1 ;
  644. X        t1 = t2 ;
  645. X    } else {
  646. X        nroots = 2 ;
  647. X    }
  648. X        
  649. X    /*
  650. X     * ensure that the points are between the two bounding planes...
  651. X     */
  652. X    
  653. X    switch(nroots) {
  654. X    case 1:
  655. X        RayPoint(ray, t1, P) ;
  656. X        d = VecDot(cd -> cone_w, P) ;
  657. X        if (d >= cd -> cone_min_d && d <= cd -> cone_max_d) {
  658. X            hit -> isect_t = t1 ;    
  659. X            hit -> isect_prim = obj ;
  660. X            hit -> isect_surf = obj -> o_surf ;
  661. X            hit -> isect_enter = 0 ;
  662. X            return(1) ;
  663. X        } else {
  664. X            return(0) ;
  665. X        }
  666. X        break ;
  667. X    case 2:
  668. X        RayPoint(ray, t1, P) ;
  669. X        d = VecDot(cd -> cone_w, P) ;
  670. X        if (d >= cd -> cone_min_d && d <= cd -> cone_max_d) {
  671. X            hit -> isect_t = t1 ;    
  672. X            hit -> isect_prim = obj ;
  673. X            hit -> isect_surf = obj -> o_surf ;
  674. X            hit -> isect_enter = 1 ;
  675. X            return(1) ;
  676. X        } else {
  677. X            RayPoint(ray, t2, P) ;
  678. X            d = VecDot(cd -> cone_w, P) ;
  679. X            if (d >= cd -> cone_min_d && d <= cd -> cone_max_d) {
  680. X                hit -> isect_t = t2 ;    
  681. X                hit -> isect_prim = obj ;
  682. X                hit -> isect_surf = obj -> o_surf ;
  683. X                hit -> isect_enter = 0 ;
  684. X                return(1) ;
  685. X            }
  686. X        }
  687. X        return(0) ;
  688. X    }
  689. X    return(0) ;
  690. X}
  691. X
  692. Xint
  693. XConeNormal(obj, hit, P, N)
  694. X Object * obj ;
  695. X Isect * hit ;
  696. X Point P, N ;
  697. X{
  698. X    Flt t ;
  699. X    Vec V ;
  700. X    ConeData * cd ;
  701. X
  702. X    cd = (ConeData *) obj -> o_data ;
  703. X
  704. X    /*
  705. X     * fill in the real normal...
  706. X     * Project the point onto the base plane.  The normal is
  707. X     * a vector from the basepoint through this point, plus the slope
  708. X     * times the cone_w vector...
  709. X     */
  710. X
  711. X    t = - (VecDot(P, cd -> cone_w) + cd -> cone_base_d) ;
  712. X    VecAddS(t, cd -> cone_w, P, V) ;
  713. X    VecSub(V, cd -> cone_base, N) ;
  714. X    VecNormalize(N) ;
  715. X    VecAddS(- cd -> cone_slope, cd -> cone_w, N, N) ;
  716. X    VecNormalize(N) ;
  717. X}
  718. X
  719. XObject *
  720. XMakeCone(basepoint, baseradius, apexpoint, apexradius) 
  721. X Vec basepoint, apexpoint ;
  722. X Flt baseradius, apexradius ;
  723. X{
  724. X    Object * obj ;
  725. X    ConeData * cd ;
  726. X    Flt dmin, dmax, d , ftmp;
  727. X    Vec tmp ;
  728. X    int i ;
  729. X
  730. X    obj = (Object *) malloc (sizeof (Object)) ;
  731. X    obj -> o_type = T_CONE ;
  732. X    obj -> o_procs = & ConeProcs ;
  733. X    obj -> o_surf = CurrentSurface ;
  734. X
  735. X    cd = (ConeData *) malloc (sizeof(ConeData)) ;
  736. X
  737. X    VecCopy(basepoint, cd -> cone_base) ;
  738. X    VecCopy(apexpoint, cd -> cone_apex) ;
  739. X
  740. X    cd -> cone_base_radius = baseradius ;
  741. X    cd -> cone_apex_radius = apexradius ;
  742. X
  743. X
  744. X    VecSub(apexpoint, basepoint, cd -> cone_w) ;
  745. X    cd -> cone_height = VecNormalize(cd -> cone_w) ;
  746. X    cd -> cone_slope =  (cd -> cone_apex_radius - cd -> cone_base_radius) /
  747. X                (cd -> cone_height) ;
  748. X    cd -> cone_base_d = - VecDot(basepoint, cd -> cone_w) ;
  749. X
  750. X    MakeVector(0.0, 0.0, 1.0, tmp) ;
  751. X
  752. X    if (1.0 - fabs(VecDot(tmp, cd ->  cone_w)) < rayeps) {
  753. X        MakeVector(0.0, 1.0, 0.0, tmp) ;
  754. X    }
  755. X
  756. X    /* find two axes which are at right angles to cone_w
  757. X     */
  758. X
  759. X    VecCross(cd -> cone_w, tmp, cd -> cone_u) ;
  760. X    VecCross(cd -> cone_u, cd -> cone_w, cd -> cone_v) ;
  761. X
  762. X    cd -> cone_min_d = VecDot(cd -> cone_w, cd -> cone_base) ;
  763. X    cd -> cone_max_d = VecDot(cd -> cone_w, cd -> cone_apex) ;
  764. X
  765. X    if (cd -> cone_max_d < cd -> cone_min_d) {
  766. X        ftmp = cd -> cone_max_d ;
  767. X        cd -> cone_max_d = cd -> cone_min_d ;
  768. X        cd -> cone_min_d = ftmp ;
  769. X    }
  770. X
  771. X    obj -> o_data = (void *) cd ;
  772. X
  773. X    for (i = 0 ; i < NSLABS ; i ++) {
  774. X        dmin = HUGE ;
  775. X        dmax = -HUGE ;
  776. X        d = VecDot(basepoint, Slab[i]) - baseradius ;
  777. X        if (d < dmin) dmin = d ; if (d > dmax) dmax = d ;
  778. X        d = VecDot(basepoint, Slab[i]) + baseradius ;
  779. X        if (d < dmin) dmin = d ; if (d > dmax) dmax = d ;
  780. X        d = VecDot(apexpoint, Slab[i]) - apexradius ;
  781. X        if (d < dmin) dmin = d ; if (d > dmax) dmax = d ;
  782. X        d = VecDot(apexpoint, Slab[i]) + apexradius ;
  783. X        if (d < dmin) dmin = d ; if (d > dmax) dmax = d ;
  784. X
  785. X        obj -> o_dmin[i] = dmin ;
  786. X        obj -> o_dmax[i] = dmax ;
  787. X    }
  788. X
  789. X    return(obj) ;
  790. X}
  791. END_OF_FILE
  792. if test 6358 -ne `wc -c <'cone.c'`; then
  793.     echo shar: \"'cone.c'\" unpacked with wrong size!
  794. fi
  795. # end of 'cone.c'
  796. fi
  797. if test -f 'intersect.c' -a "${1}" != "-c" ; then 
  798.   echo shar: Will not clobber existing file \"'intersect.c'\"
  799. else
  800. echo shar: Extracting \"'intersect.c'\" \(4725 characters\)
  801. sed "s/^X//" >'intersect.c' <<'END_OF_FILE'
  802. X/***********************************************************************
  803. X * $Author: markv $
  804. X * $Revision: 1.2 $
  805. X * $Date: 88/09/12 12:54:48 $
  806. X * $Log:    intersect.c,v $
  807. X * Revision 1.2  88/09/12  12:54:48  markv
  808. X * Added early cutoffs, as suggested by Haines in the RT-News, and 
  809. X * independantly discovered by myself during our correspondence.
  810. X * Now, Intersect takes a max distance, and will not search for 
  811. X * any intersections beyond the max distance.
  812. X * Also, to enable shadow caching, Shadow now returns the object
  813. X * that it found on its way to the light source.
  814. X * 
  815. X * These optimizations speed up shadow testing considerably, and 
  816. X * normal rays some small amount.
  817. X * 
  818. X * Revision 1.1  88/09/11  11:00:40  markv
  819. X * Initial revision
  820. X * 
  821. X ***********************************************************************/
  822. X
  823. X#include <stdio.h>
  824. X#include <math.h>
  825. X#include <assert.h>
  826. X#include "defs.h"
  827. X#include "extern.h"
  828. X
  829. X/*
  830. X * intersect.c
  831. X * Much nicer now, uses the nifty priority queue search
  832. X * as suggested by Kajiya...
  833. X */
  834. X
  835. XFlt        num[NSLABS] ;
  836. XFlt        den[NSLABS] ;
  837. X
  838. X/***********************************************************************
  839. X * CheckAndEnqueue(obj, maxdist)
  840. X * Check the current ray (as paramaterized with the num and den 
  841. X * arrays above) against the bounding volume of obj.
  842. X * If we intersect the bounding volume, then insert it into the 
  843. X * priority queue.
  844. X *
  845. X * Note: should be broken into two separate procedures...
  846. X ***********************************************************************/
  847. X
  848. XINLINE
  849. XCheckAndEnqueue(obj, maxdist)
  850. X Object * obj ;
  851. X Flt maxdist ;
  852. X{
  853. X    int i ;
  854. X    Flt tmp ;
  855. X    Flt tmin, tmax ;
  856. X    Flt dmin = -HUGE ;
  857. X    Flt dmax = maxdist ;
  858. X
  859. X    nChecked ++ ;
  860. X
  861. X    for (i = 0 ; i < NSLABS ; i ++) {
  862. X
  863. X        /* enters the slab here...    */
  864. X        tmin = (obj -> o_dmin[i] - num[i]) * den[i] ;
  865. X        /* and exits here...        */
  866. X        tmax = (obj -> o_dmax[i] - num[i]) * den[i] ;
  867. X
  868. X        /* but we may have to swap...    */
  869. X        if (tmin > tmax) {
  870. X            tmp = tmin ; tmin = tmax ; tmax = tmp ;
  871. X        }
  872. X
  873. X        /* if exited closer than we thought, update     */
  874. X        if (tmax < dmax)
  875. X            dmax = tmax ;
  876. X        /* if entered farther than we thought, update     */
  877. X        if (tmin > dmin)
  878. X            dmin = tmin ;
  879. X
  880. X        if (dmin > dmax || dmax < rayeps)
  881. X            return ;
  882. X    }
  883. X    PriorityQueueInsert(dmin, obj) ;
  884. X    nEnqueued ++ ;
  885. X}
  886. X
  887. X/***********************************************************************
  888. X * Intersect(ray, hit, maxdist)
  889. X * 
  890. X * Returns true if we hit something in the root model closer than maxdist.  
  891. X * Returns the closest hit in the "hit" buffer.
  892. X ***********************************************************************/
  893. X
  894. XIntersect(ray, hit, maxdist)
  895. X Ray * ray ;
  896. X Isect * hit ;
  897. X Flt maxdist ;
  898. X{
  899. X    Isect        nhit ;
  900. X    int        i ;
  901. X    Flt        min_dist = maxdist ;
  902. X    Object *    cobj ;
  903. X    Object *     pobj = NULL ;
  904. X    CompositeData     * cdp ;
  905. X    Flt        key ;
  906. X
  907. X    /* If the object is simple, then return the hit that it gives you */
  908. X
  909. X    if (Root -> o_type != T_COMPOSITE)
  910. X        return (Root -> o_procs -> intersect) (Root, ray, hit) ;
  911. X
  912. X    for (i = 0 ; i < NSLABS ; i ++) {
  913. X        num[i] = VecDot(ray -> P, Slab[i]) ;
  914. X        den[i] = 1.0 / VecDot(ray -> D, Slab[i]) ;
  915. X    }
  916. X
  917. X    /* start with an empty priority queue */
  918. X    PriorityQueueNull() ;
  919. X
  920. X    CheckAndEnqueue(Root, maxdist) ;
  921. X
  922. X    for (;;) {
  923. X
  924. X        if (PriorityQueueEmpty())
  925. X            break ;
  926. X
  927. X        PriorityQueueDelete(&key, &cobj) ;
  928. X
  929. X        if (key > min_dist) {
  930. X
  931. X            /*
  932. X             * we have already found a primitive
  933. X             * that was closer, we need look no further...
  934. X             */
  935. X             break ;
  936. X
  937. X        } else if (cobj -> o_type == T_COMPOSITE) {
  938. X            /* 
  939. X             * if it is in the queue, it got hit.
  940. X             * check each of its children to see if their
  941. X             * bounding volumes get hit.
  942. X             * if so, then push them into the priority
  943. X             * queue...
  944. X             */
  945. X            
  946. X            cdp = (CompositeData *) cobj -> o_data ;
  947. X
  948. X            for (i = 0 ; i < cdp -> c_size ; i ++ ) {
  949. X                CheckAndEnqueue(cdp -> c_object[i], maxdist) ;
  950. X            }
  951. X
  952. X        } else {
  953. X
  954. X            /*
  955. X             * we have a primitive 
  956. X             * intersect with the primitive, and possibly
  957. X             * update the nearest hit if it is indeed closer
  958. X             * than the one we currently have...
  959. X             */
  960. X            
  961. X            if ((cobj -> o_procs -> intersect) (cobj, ray, &nhit)) {
  962. X                if (nhit.isect_t < min_dist) {
  963. X                    pobj = cobj ;
  964. X                    *hit = nhit ;
  965. X                    min_dist = nhit.isect_t ;
  966. X                }
  967. X            }
  968. X        }
  969. X    }
  970. X
  971. X    if (pobj)
  972. X        return 1 ;
  973. X    else
  974. X        return 0 ;
  975. X}
  976. X
  977. X/***********************************************************************
  978. X * Shadow(ray, hit, tmax)
  979. X * 
  980. X * Returns true if we are unshadowed.  Returns the primitive in 
  981. X * the "hit" buffer.
  982. X *
  983. X * Note: the return value of this procedure is a bit strange, as well 
  984. X * as the name.  Should probably be changed.
  985. X ***********************************************************************/
  986. X
  987. Xint 
  988. XShadow(ray, hit, tmax) 
  989. X Ray * ray ;
  990. X Isect * hit ;
  991. X Flt tmax ;
  992. X{
  993. X    if (Intersect(ray, hit, tmax))
  994. X        return 0 ;
  995. X    else
  996. X        return 1 ;
  997. X}
  998. END_OF_FILE
  999. if test 4725 -ne `wc -c <'intersect.c'`; then
  1000.     echo shar: \"'intersect.c'\" unpacked with wrong size!
  1001. fi
  1002. # end of 'intersect.c'
  1003. fi
  1004. if test -f 'poly.c' -a "${1}" != "-c" ; then 
  1005.   echo shar: Will not clobber existing file \"'poly.c'\"
  1006. else
  1007. echo shar: Extracting \"'poly.c'\" \(5428 characters\)
  1008. sed "s/^X//" >'poly.c' <<'END_OF_FILE'
  1009. X/***********************************************************************
  1010. X * $Author: markv $
  1011. X * $Revision: 1.2 $
  1012. X * $Date: 88/10/04 14:32:25 $
  1013. X * $Log:    poly.c,v $
  1014. X * Revision 1.2  88/10/04  14:32:25  markv
  1015. X * Renamed p1 and p2 to be poly_u and poly_v, which are better names.
  1016. X * Also may have solved problems having to do with floating point roundoff
  1017. X * when planes are parallel to bounding slabs.
  1018. X * 
  1019. X * Revision 1.1  88/09/11  11:00:42  markv
  1020. X * Initial revision
  1021. X * 
  1022. X ***********************************************************************/
  1023. X
  1024. X#include <stdio.h>
  1025. X#include <math.h>
  1026. X#include "defs.h"
  1027. X#include "extern.h"
  1028. X
  1029. Xtypedef struct t_polydata {
  1030. X    int     poly_npoints ;
  1031. X    Vec    * poly_point ;
  1032. X    Vec    poly_normal ;
  1033. X    Flt     poly_d ;
  1034. X    Flt    poly_u, poly_v ;
  1035. X} PolyData ;
  1036. X
  1037. Xint PolyPrint ();
  1038. Xint PolyIntersect ();
  1039. Xint PolyNormal ();
  1040. X
  1041. XObjectProcs PolyProcs = {
  1042. X    PolyPrint,
  1043. X    PolyIntersect,
  1044. X    PolyNormal,
  1045. X} ;
  1046. X
  1047. Xint 
  1048. XPolyPrint(obj)
  1049. X Object * obj ;
  1050. X{
  1051. X    int i ;
  1052. X    PolyData * pd ;
  1053. X
  1054. X    pd = (PolyData *) obj -> o_data ;
  1055. X    printf("p %d\n", pd -> poly_npoints) ;
  1056. X    for (i = 0 ; i < pd -> poly_npoints ; i++) {
  1057. X        printf("%g %g %g\n", pd -> poly_point[i][0],
  1058. X                    pd -> poly_point[i][1],
  1059. X                    pd -> poly_point[i][2]) ;
  1060. X    }
  1061. X}
  1062. X
  1063. X/***********************************************************************
  1064. X * PolyIntersect(obj, ray, hit)
  1065. X * 
  1066. X * returns 1 if we hit the polygon, with the hit information in hit.
  1067. X * Uses a version of Jordan's theorem to determine whether the point 
  1068. X * is inside the polygon.
  1069. X * The variable "crosses" will count the number of times that we
  1070. X * cross the boundary of the curve.  If it is odd, we are inside.
  1071. X *
  1072. X ***********************************************************************/
  1073. X
  1074. XPolyIntersect(obj, ray, hit)
  1075. X Object * obj ;
  1076. X Ray * ray ;
  1077. X Isect * hit ;
  1078. X{
  1079. X    Flt n,d,t,m,b ;
  1080. X    Point V ;
  1081. X    int i, j, crosses = 0 ;
  1082. X    int qi, qj ;
  1083. X    int ri, rj ;
  1084. X    int u, v ;
  1085. X    PolyData * pd ;
  1086. X
  1087. X    pd = (PolyData *) obj -> o_data ;
  1088. X    n = VecDot(ray -> P, pd -> poly_normal) + pd -> poly_d ;
  1089. X    d = VecDot(ray -> D, pd -> poly_normal) ;
  1090. X
  1091. X    if ((Flt) fabs(d) < rayeps) {
  1092. X        return (0) ;
  1093. X    }
  1094. X
  1095. X    t = - n / d ;
  1096. X    if (t < rayeps)
  1097. X        return 0 ;
  1098. X
  1099. X    RayPoint(ray,t,V);
  1100. X
  1101. X    u = pd -> poly_u ;
  1102. X    v = pd -> poly_v ;
  1103. X
  1104. X    for (i = 0 ; i < pd -> poly_npoints ; i++) {
  1105. X
  1106. X        j = (i + 1) % pd -> poly_npoints ;
  1107. X
  1108. X        qi = 0 ; qj = 0 ;
  1109. X        ri = 0 ; rj = 0 ;
  1110. X
  1111. X        if (pd -> poly_point[i][v] == pd -> poly_point[j][v])
  1112. X            continue ;        /*ignore horizontal lines */
  1113. X
  1114. X        /*
  1115. X         * If we are both above, or both below the intersection point,
  1116. X         * go onto the next one.
  1117. X         */
  1118. X
  1119. X        if (pd -> poly_point[i][v] < V[v])
  1120. X            qi = 1 ;
  1121. X        if (pd -> poly_point[j][v] < V[v])
  1122. X            qj = 1 ;
  1123. X        if (qi == qj)
  1124. X            continue ;
  1125. X
  1126. X        /*
  1127. X         * We know one end point was above, and one was below.
  1128. X         * If theyare both to the left, then we crossed the line 
  1129. X         * to negative infinity, and we continue.
  1130. X         */
  1131. X        if (pd -> poly_point[i][u] < V[u])
  1132. X            ri = 1 ;
  1133. X         if (pd -> poly_point[j][u] < V[u])
  1134. X            rj = 1 ;
  1135. X
  1136. X        if (ri & rj) {
  1137. X            crosses ++ ;
  1138. X            continue ;
  1139. X        }
  1140. X
  1141. X        /*
  1142. X         * Otherwise, if we are both to the right, 
  1143. X         * we can continue without a cross.
  1144. X         */
  1145. X
  1146. X        if ((ri|rj) == 0)
  1147. X            continue ;
  1148. X
  1149. X
  1150. X        /* 
  1151. X         * more difficult acceptance...
  1152. X         * We have a line segment which occurs with endpoints
  1153. X         * in diagonally opposite quadrants.  We must solve
  1154. X         * for the intersection, ie, where v = 0.
  1155. X         */
  1156. X
  1157. X        m = (pd -> poly_point[j][v] - pd -> poly_point[i][v]) /
  1158. X            (pd -> poly_point[j][u] - pd -> poly_point[i][u]) ;
  1159. X        
  1160. X        b = (pd -> poly_point[j][v] - V[v]) - 
  1161. X            m * (pd -> poly_point[j][u] - V[u]);
  1162. X
  1163. X        if ((-b/m) < rayeps)
  1164. X            crosses ++ ;
  1165. X    }
  1166. X
  1167. X    if (crosses & 1) {
  1168. X        hit -> isect_t = t ;
  1169. X        hit -> isect_surf = obj -> o_surf ;
  1170. X        hit -> isect_prim = obj ;
  1171. X        hit -> isect_enter = 0 ;
  1172. X        return(1);
  1173. X    } else {
  1174. X        return(0);
  1175. X    }
  1176. X}
  1177. X
  1178. XPolyNormal(obj, hit, P, N)
  1179. X Object * obj ;
  1180. X Isect * hit ;
  1181. X Point P, N ;
  1182. X{
  1183. X
  1184. X    PolyData * pd ;
  1185. X    pd = (PolyData *) obj -> o_data ;
  1186. X    VecCopy(pd -> poly_normal, N);
  1187. X}
  1188. X
  1189. XObject *
  1190. XMakePoly(npoints, points)
  1191. X int npoints ;
  1192. X Vec * points ;
  1193. X{
  1194. X    Object * obj ;
  1195. X    PolyData * pd ;
  1196. X    Vec P1, P2 ;
  1197. X    Flt d, dmax, dmin ;
  1198. X    int i, j ;
  1199. X
  1200. X    obj = (Object *) malloc (sizeof(Object)) ;
  1201. X    obj -> o_type = T_POLY ;
  1202. X    obj -> o_procs = & PolyProcs ;
  1203. X    obj -> o_surf = CurrentSurface ;
  1204. X
  1205. X    pd = (PolyData *) malloc (sizeof(PolyData)) ;
  1206. X    pd -> poly_npoints = npoints ;
  1207. X    pd -> poly_point = points ;
  1208. X
  1209. X    /*
  1210. X     * calculate the normal by giving various cross products...
  1211. X     */
  1212. X    
  1213. X    VecSub(pd -> poly_point[0], pd -> poly_point[1], P1) ;
  1214. X    VecSub(pd -> poly_point[2], pd -> poly_point[1], P2) ;
  1215. X
  1216. X    VecCross(P1, P2, pd -> poly_normal) ;
  1217. X    VecNormalize(pd -> poly_normal) ;
  1218. X
  1219. X    if (fabs(pd -> poly_normal[0]) > fabs(pd -> poly_normal[1])
  1220. X        && fabs(pd -> poly_normal[0]) > fabs(pd -> poly_normal[2])) {
  1221. X        pd -> poly_u = 1 ;
  1222. X        pd -> poly_v = 2 ;
  1223. X    } else if (fabs(pd -> poly_normal[1]) > fabs(pd -> poly_normal[0]) 
  1224. X        && fabs(pd -> poly_normal[1]) > fabs(pd -> poly_normal[2])) {
  1225. X        pd -> poly_u = 0 ;
  1226. X        pd -> poly_v = 2 ;
  1227. X    } else {
  1228. X        pd -> poly_u = 0 ;
  1229. X        pd -> poly_v = 1 ;
  1230. X    }
  1231. X
  1232. X    pd -> poly_d = - VecDot(pd -> poly_normal, pd -> poly_point[0]) ;
  1233. X
  1234. X    obj -> o_data = (void *) pd ;
  1235. X
  1236. X    /*
  1237. X     * now, calculate the values of 
  1238. X     * the dmin and dmax 'es for the globally defined slabs...
  1239. X     */
  1240. X    
  1241. X    for (i = 0 ; i < NSLABS ; i ++) {
  1242. X        dmin = HUGE ;
  1243. X        dmax = - HUGE ;
  1244. X
  1245. X        for (j = 0 ; j < pd -> poly_npoints ; j ++) {
  1246. X            d = VecDot(Slab[i], pd -> poly_point[j]) ;
  1247. X            if (d < dmin) dmin = d ;
  1248. X            if (d > dmax) dmax = d ;
  1249. X        }
  1250. X        obj -> o_dmin[i] = dmin - rayeps ;
  1251. X        obj -> o_dmax[i] = dmax + rayeps ;
  1252. X    }
  1253. X    return(obj) ;
  1254. X}
  1255. END_OF_FILE
  1256. if test 5428 -ne `wc -c <'poly.c'`; then
  1257.     echo shar: \"'poly.c'\" unpacked with wrong size!
  1258. fi
  1259. # end of 'poly.c'
  1260. fi
  1261. if test -f 'screen.c' -a "${1}" != "-c" ; then 
  1262.   echo shar: Will not clobber existing file \"'screen.c'\"
  1263. else
  1264. echo shar: Extracting \"'screen.c'\" \(5370 characters\)
  1265. sed "s/^X//" >'screen.c' <<'END_OF_FILE'
  1266. X/***********************************************************************
  1267. X * $Author: markv $
  1268. X * $Revision: 1.3 $
  1269. X * $Date: 88/10/04 14:33:37 $
  1270. X * $Log:    screen.c,v $
  1271. X * Revision 1.3  88/10/04  14:33:37  markv
  1272. X * Revamped most of this code.  Added code for no antialiasing, jittering, 
  1273. X * filtered, and statistically optimized antialiasing.   Statistically
  1274. X * optimized anti aliasing still doesn't work properly.
  1275. X * 
  1276. X * Revision 1.2  88/09/12  10:00:54  markv
  1277. X * Fixed lingering (possible) bug with curbuf memory allocation.  
  1278. X * Size returned by malloc is now correct in both spots.
  1279. X * Thanks tim@ben.
  1280. X * 
  1281. X * Revision 1.1  88/09/11  11:00:43  markv
  1282. X * Initial revision
  1283. X * 
  1284. X ***********************************************************************/
  1285. X
  1286. X#include <stdio.h>
  1287. X#include <math.h>
  1288. X#include <assert.h>
  1289. X#include "defs.h"
  1290. X#include "pic.h"
  1291. X#include "extern.h"
  1292. X
  1293. XScreen(view, picfile, xres, yres)
  1294. X Viewpoint *view ;
  1295. X char *picfile ;
  1296. X int xres, yres ;
  1297. X{
  1298. X    Pixel     *buf, *oldbuf, *curbuf, *tmp ;
  1299. X    Pic      *pic, *PicOpen();
  1300. X    Point     viewvec, leftvec, upvec ;
  1301. X    int red, green, blue ;
  1302. X    Point    dir ;
  1303. X    Ray     ray ;
  1304. X    Color    color;
  1305. X    Flt    frustrumwidth ;
  1306. X    int     x, y ;
  1307. X
  1308. X    pic = PicOpen(picfile, xres, yres) ;
  1309. X
  1310. X    /*
  1311. X     * calculate the "up" vector...
  1312. X     */
  1313. X
  1314. X    VecCopy(view -> view_up, upvec) ;
  1315. X    VecNormalize(upvec) ;
  1316. X
  1317. X    /*
  1318. X     * and the "view" vector....
  1319. X     */
  1320. X
  1321. X    VecSub(view -> view_at, view-> view_from, viewvec);
  1322. X    VecNormalize(viewvec);
  1323. X
  1324. X    /*
  1325. X     * And the "left" vector...
  1326. X     */
  1327. X
  1328. X    VecCross(view -> view_up, viewvec, leftvec);
  1329. X    VecNormalize(leftvec);
  1330. X
  1331. X    /*
  1332. X     * Calculate the width of the view frustrum in world coordinates.
  1333. X     * and then scale the left and up vectors appropriately.
  1334. X     */
  1335. X
  1336. X    frustrumwidth = (view -> view_dist) * ((Flt) tan(view -> view_angle)) ;
  1337. X    VecScale(-frustrumwidth, upvec) ;
  1338. X    VecScale(-frustrumwidth, leftvec) ;
  1339. X
  1340. X    if (filtflag) {
  1341. X        FilterScan( pic,
  1342. X            view -> view_from,
  1343. X            viewvec,
  1344. X            upvec,
  1345. X            leftvec, 
  1346. X            xres, yres) ;
  1347. X    } else if (jitterflag) {
  1348. X        JitterScan(   pic,
  1349. X            view -> view_from,
  1350. X            viewvec,
  1351. X            upvec,
  1352. X            leftvec, 
  1353. X             xres, yres) ;
  1354. X    } else {
  1355. X        Scan(   pic,
  1356. X            view -> view_from,
  1357. X            viewvec,
  1358. X            upvec,
  1359. X            leftvec, 
  1360. X             xres, yres) ;
  1361. X    }
  1362. X
  1363. X    PicClose(pic) ;
  1364. X}
  1365. X
  1366. X
  1367. X
  1368. X
  1369. XScan(pic, eye, viewvec, upvec, leftvec, xres, yres)
  1370. X Pic * pic ;
  1371. X Vec eye ;
  1372. X Vec viewvec ;
  1373. X Vec upvec ;
  1374. X Vec leftvec ;
  1375. X int xres, yres ;
  1376. X{
  1377. X    Ray ray ;
  1378. X    int x, y ;
  1379. X    Flt xlen, ylen ;
  1380. X    Color color ;
  1381. X
  1382. X    /*
  1383. X     * Generate the image...
  1384. X     */
  1385. X
  1386. X    VecCopy(eye, ray.P) ;
  1387. X    for (y = 0 ; y < yres ; y++) {
  1388. X        for (x = 0 ; x < xres ; x++) {
  1389. X            xlen = ((Flt) (2 * x) / (Flt) xres) - 1.0 ;
  1390. X            ylen = ((Flt) (2 * y) / (Flt) yres) - 1.0 ;
  1391. X            VecComb(xlen, leftvec, ylen, upvec, ray.D) ;
  1392. X            VecAdd(ray.D, viewvec, ray.D) ;
  1393. X            VecNormalize(ray.D);
  1394. X            Trace(0, 1.0, &ray, color);
  1395. X            color[0] = min(1.0, color[0]) ;
  1396. X            color[1] = min(1.0, color[1]) ;
  1397. X            color[2] = min(1.0, color[2]) ;
  1398. X            PicWritePixel(pic, color) ;
  1399. X        }
  1400. X        if (tickflag)
  1401. X            fprintf(stderr, ".") ;
  1402. X    }
  1403. X    if (tickflag)
  1404. X        fprintf(stderr, "\n") ;
  1405. X}
  1406. X
  1407. XFilterScan(pic, eye, viewvec, upvec, leftvec, xres, yres)
  1408. X Pic * pic ;
  1409. X Vec eye ;
  1410. X Vec viewvec ;
  1411. X Vec upvec ;
  1412. X Vec leftvec ;
  1413. X int xres, yres ;
  1414. X{
  1415. X    Ray ray ;
  1416. X    int x, y, i ;
  1417. X    Flt xlen, ylen ;
  1418. X    Color color ;
  1419. X    Color * nbuf, *obuf, *tmp  ;    /* new and old buffers, resp.    */
  1420. X    Color avg ;
  1421. X
  1422. X    /*
  1423. X     * allocate enough memory for the filter buffer
  1424. X     */
  1425. X
  1426. X    nbuf = (Color *) calloc ((xres + 1), sizeof(Color)) ;
  1427. X    obuf = NULL ;
  1428. X
  1429. X    VecCopy(eye, ray.P) ;
  1430. X
  1431. X    for (y = 0 ; y <= yres ; y++) {
  1432. X        for (x = 0 ; x <= xres ; x++) {
  1433. X            xlen = ((Flt) (2 * x) / (Flt) xres) - 1.0 ;
  1434. X            ylen = ((Flt) (2 * y) / (Flt) yres) - 1.0 ;
  1435. X            VecComb(xlen, leftvec, ylen, upvec, ray.D) ;
  1436. X            VecAdd(ray.D, viewvec, ray.D) ;
  1437. X            VecNormalize(ray.D);
  1438. X            Trace(0, 1.0, &ray, color);
  1439. X            color[0] = min(1.0, color[0]) ;
  1440. X            color[1] = min(1.0, color[1]) ;
  1441. X            color[2] = min(1.0, color[2]) ;
  1442. X            VecCopy(color, nbuf[x]) ;
  1443. X        }
  1444. X        if (obuf) {
  1445. X            for (i = 0 ; i < xres ; i++) {
  1446. X                VecAdd(nbuf[i], nbuf[i+1], avg) ;
  1447. X                VecAdd(obuf[i], avg, avg) ;
  1448. X                VecAdd(obuf[i+1], avg, avg) ;
  1449. X                VecScale(0.25, avg) ;
  1450. X                PicWritePixel(pic, avg) ;
  1451. X            }
  1452. X            tmp = obuf ;
  1453. X            obuf = nbuf ;
  1454. X            nbuf = tmp ;
  1455. X        } else {
  1456. X            /* 
  1457. X             * first scan line, set it up wierdly...
  1458. X             */
  1459. X            obuf = nbuf ;
  1460. X            nbuf = (Color *) calloc ((xres + 1), sizeof(Color)) ;
  1461. X        }
  1462. X        if (tickflag)
  1463. X            fprintf(stderr, ".") ;
  1464. X    }
  1465. X    if (tickflag)
  1466. X        fprintf(stderr, "\n") ;
  1467. X}
  1468. X
  1469. XJitterScan(pic, eye, viewvec, upvec, leftvec, xres, yres)
  1470. X Pic * pic ;
  1471. X Vec eye ; 
  1472. X Vec viewvec ; 
  1473. X Vec upvec ; 
  1474. X Vec leftvec ;
  1475. X int xres, yres ;
  1476. X{
  1477. X    Ray ray ;
  1478. X    int x, y, i ;
  1479. X    Flt xlen, ylen ;
  1480. X    Flt xwidth, ywidth ;
  1481. X    Color color, avg ;
  1482. X
  1483. X    /*
  1484. X     * Generate the image...
  1485. X     */
  1486. X
  1487. X    xwidth = 2.0 / (Flt) xres ;
  1488. X    ywidth = 2.0 / (Flt) xres ;
  1489. X
  1490. X    VecCopy(eye, ray.P) ;
  1491. X    for (y = 0 ; y < yres ; y++) {
  1492. X        ylen = ((Flt) (2 * y) / (Flt) yres) - 1.0 ;
  1493. X        for (x = 0 ; x < xres ; x++) {
  1494. X            xlen = ((Flt) (2 * x) / (Flt) xres) - 1.0 ;
  1495. X
  1496. X            avg[0] = avg[1] = avg[2] = 0.0 ;
  1497. X
  1498. X            for (i = 0 ; i < maxsamples ; i++) {
  1499. X                VecComb(xlen + rnd() * xwidth, leftvec, 
  1500. X                    ylen + rnd() * ywidth, upvec, ray.D) ;
  1501. X                VecAdd(ray.D, viewvec, ray.D) ;
  1502. X                VecNormalize(ray.D);
  1503. X                Trace(0, 1.0, &ray, color);
  1504. X                VecAdd(color, avg, avg) ;
  1505. X            }
  1506. X
  1507. X            VecScale(1.0 / (Flt) maxsamples, avg) ;
  1508. X            avg[0] = min(1.0, avg[0]) ;
  1509. X            avg[1] = min(1.0, avg[1]) ;
  1510. X            avg[2] = min(1.0, avg[2]) ;
  1511. X            PicWritePixel(pic, avg) ;
  1512. X        }
  1513. X        if (tickflag)
  1514. X            fprintf(stderr, ".") ;
  1515. X    }
  1516. X    if (tickflag)
  1517. X        fprintf(stderr, "\n") ;
  1518. X}
  1519. END_OF_FILE
  1520. if test 5370 -ne `wc -c <'screen.c'`; then
  1521.     echo shar: \"'screen.c'\" unpacked with wrong size!
  1522. fi
  1523. # end of 'screen.c'
  1524. fi
  1525. if test -f 'shade.c' -a "${1}" != "-c" ; then 
  1526.   echo shar: Will not clobber existing file \"'shade.c'\"
  1527. else
  1528. echo shar: Extracting \"'shade.c'\" \(4575 characters\)
  1529. sed "s/^X//" >'shade.c' <<'END_OF_FILE'
  1530. X/***********************************************************************
  1531. X * $Author: markv $
  1532. X * $Revision: 1.2 $
  1533. X * $Date: 88/09/12 12:58:01 $
  1534. X * $Log:    shade.c,v $
  1535. X * Revision 1.2  88/09/12  12:58:01  markv
  1536. X * Added specular reflections and shadow caching.  When a object
  1537. X * is found between a light source and the current point we are trying
  1538. X * to shade, that object is cached (indexed by recursion level) in the
  1539. X * lightsource.  Next time we test a shadow against the light source,
  1540. X * we test this object first.  If it is between us and the light source,
  1541. X * we can correctly shadow the object without calling Intersect().
  1542. X * 
  1543. X * Note: specular highlights call the unix pow() function, which seems
  1544. X * to be REALLY expensive.  Optimizations could be made here.
  1545. X * 
  1546. X * Revision 1.1  88/09/11  11:00:44  markv
  1547. X * Initial revision
  1548. X *  
  1549. X ***********************************************************************/
  1550. X
  1551. X#include <stdio.h>
  1552. X#include <math.h>
  1553. X#include "defs.h"
  1554. X#include "extern.h"
  1555. X
  1556. X/***********************************************************************
  1557. X * Shade(level, weight, P, N, I, hit, col)
  1558. X *
  1559. X * Wow! Too many parameters!
  1560. X *
  1561. X * Shade is the driving force of the raytracer.  It calculates the actual
  1562. X * luminance at point P, given N, the normal, and I the incident vector.
  1563. X * We also pass in the hit, because the normal generating code might need it.
  1564. X * The color is returned in col.
  1565. X ***********************************************************************/
  1566. X
  1567. XShade(level, weight, P, N, I, hit, col) 
  1568. X int level ;
  1569. X Flt weight ;
  1570. X Vec P, N, I ;
  1571. X Isect * hit;
  1572. X Color col ;
  1573. X{
  1574. X    Ray     tray ;
  1575. X    Color     tcol ;
  1576. X    Vec     L, H, R ;
  1577. X    Flt     t ;
  1578. X    Flt    diff ;
  1579. X    Flt     spec ;
  1580. X#ifdef SHADOW_CACHING
  1581. X    Object    * cached ;
  1582. X#endif /* SHADOW_CACHING */
  1583. X    Isect nhit ;
  1584. X    Surface    * surf ;
  1585. X    int l ;
  1586. X
  1587. X    col[0] = col[1] = col[2] = 0.0 ;
  1588. X    surf = hit -> isect_prim -> o_surf ;
  1589. X
  1590. X    for (l = 0; l < nLights; l++) {
  1591. X        VecSub(Lights[l].light_pos, P, L);
  1592. X        if (VecDot(N,L) >= 0.0) {
  1593. X            t = VecNormalize(L);
  1594. X            VecCopy(P, tray.P);
  1595. X            VecCopy(L, tray.D);
  1596. X            nShadows ++ ;
  1597. X#ifdef SHADOW_CACHING
  1598. X            cached = Lights[l].light_obj_cache[level] ;
  1599. X            if (cached 
  1600. X                && (cached -> o_procs -> intersect) 
  1601. X                    (cached, &tray, &nhit) 
  1602. X                && nhit.isect_t < t) {
  1603. X                /* 
  1604. X                 * we are in shadow, continue...
  1605. X                 */
  1606. X                nShadowCacheHits ++ ;
  1607. X                continue ;
  1608. X            }
  1609. X#endif /* SHADOW_CACHING */
  1610. X            if (Shadow(&tray, &nhit, t)) {
  1611. X                diff = VecDot(N,L) * surf -> surf_kd 
  1612. X                    * Lights[l].light_brightness ;
  1613. X#ifdef SHADOW_CACHING
  1614. X                Lights[l].light_obj_cache[level] = NULL ;
  1615. X#endif /* SHADOW_CACHING */
  1616. X                VecAddS(diff, surf -> surf_color, col, col) ;
  1617. X                SpecularDirection(I, N, R) ;
  1618. X                VecNormalize(R) ;
  1619. X                if (surf -> surf_shine > rayeps) {
  1620. X                    spec = VecDot(R,L) ;
  1621. X                    if (spec > rayeps) {
  1622. X                        spec = pow(spec, surf -> surf_shine ) * Lights[l].light_brightness ;
  1623. X                        col[0] += spec ; 
  1624. X                        col[1] += spec ; 
  1625. X                        col[2] += spec ;
  1626. X                    }
  1627. X                }
  1628. X            } else {
  1629. X#ifdef SHADOW_CACHING
  1630. X                Lights[l].light_obj_cache[level] = 
  1631. X                    nhit.isect_prim ;
  1632. X#endif /* SHADOW_CACHING */
  1633. X            }
  1634. X                
  1635. X        }
  1636. X    }
  1637. X
  1638. X    VecCopy(P, tray.P);
  1639. X
  1640. X    if(surf -> surf_ks * weight > minweight) {
  1641. X        nReflected ++ ;
  1642. X        SpecularDirection(I, N, tray.D);
  1643. X        VecNormalize(tray.D);
  1644. X        Trace(level + 1, surf -> surf_ks * weight, &tray, tcol);
  1645. X        VecAddS(surf -> surf_ks, tcol, col, col);
  1646. X    }
  1647. X
  1648. X    if (surf -> surf_kt * weight > minweight) { 
  1649. X        nRefracted ++ ;
  1650. X        if (hit -> isect_enter) 
  1651. X            TransmissionDirection(NULL, surf, I, N, tray.D) ;
  1652. X        else    
  1653. X            TransmissionDirection(surf, NULL, I, N, tray.D) ;
  1654. X        Trace(level + 1, surf -> surf_kt * weight, &tray, tcol) ;
  1655. X        VecAddS(surf -> surf_kt, tcol, col, col) ;
  1656. X    }
  1657. X}
  1658. X
  1659. X/***********************************************************************
  1660. X * SpecularDirection(I, N, R)
  1661. X * 
  1662. X * Given an incident vector I, and the normal N, calculate the 
  1663. X * direction of the reflected ray R.
  1664. X ***********************************************************************/
  1665. X
  1666. XSpecularDirection(I, N, R)
  1667. X Vec I, N, R;
  1668. X{
  1669. X    VecComb(1.0/fabs(VecDot(I,N)), I, 2.0, N, R);
  1670. X    VecNormalize(R);
  1671. X}
  1672. X
  1673. X/***********************************************************************
  1674. X * TransmissionDirection(m1, m2, I, N, T)
  1675. X *
  1676. X * calculates the direction of the transmitted ray
  1677. X ***********************************************************************/
  1678. X
  1679. XTransmissionDirection(m1, m2, I, N, T)
  1680. X Surface *m1, *m2;
  1681. X Vec I, N, T ;
  1682. X{
  1683. X    Flt n1, n2, eta, c1, cs2 ;
  1684. X    n1 = m1 ? m1 -> surf_ior : 1.0 ;
  1685. X    n2 = m2 ? m2 -> surf_ior : 1.0 ;
  1686. X    eta = n1/n2 ;
  1687. X
  1688. X    c1 = -VecDot(I,N);
  1689. X    cs2 = 1.0 - eta * eta*(1.0 - c1*c1);
  1690. X    if (cs2 < 0.0)
  1691. X        return 0;
  1692. X    VecComb(eta, I, eta*c1-sqrt(cs2), N, T);
  1693. X    return(1);
  1694. X}
  1695. END_OF_FILE
  1696. if test 4575 -ne `wc -c <'shade.c'`; then
  1697.     echo shar: \"'shade.c'\" unpacked with wrong size!
  1698. fi
  1699. # end of 'shade.c'
  1700. fi
  1701. if test -f 'tri.c' -a "${1}" != "-c" ; then 
  1702.   echo shar: Will not clobber existing file \"'tri.c'\"
  1703. else
  1704. echo shar: Extracting \"'tri.c'\" \(6022 characters\)
  1705. sed "s/^X//" >'tri.c' <<'END_OF_FILE'
  1706. X/***********************************************************************
  1707. X * $Author: markv $
  1708. X * $Revision: 1.1 $
  1709. X * $Date: 88/11/12 16:23:13 $
  1710. X * $Log:    tri.c,v $
  1711. X * Revision 1.1  88/11/12  16:23:13  markv
  1712. X * Initial revision
  1713. X * 
  1714. X * TRIs are triangular patches with normals defined at the vertices.
  1715. X * When an intersection is found, it interpolates the normal to the 
  1716. X * surface at that point.
  1717. X *
  1718. X * Algorithm is due to Jeff Arenburg, (arenburg@trwrb.uucp) and was 
  1719. X * posted to USENET.
  1720. X *
  1721. X * Basically, for each triangle we calculate an inverse transformation
  1722. X * matrix, and use it to determine the point of intersection in the plane
  1723. X * of the triangle relative to the "base point" of the triangle.  We then
  1724. X * figure its coordinates relative to that base point.  These base points
  1725. X * are used to find the barycentric coordinates, and then in normal 
  1726. X * interpolation...
  1727. X ***********************************************************************/
  1728. X
  1729. X#include <stdio.h>
  1730. X#include <math.h>
  1731. X#include "defs.h"
  1732. X#include "extern.h"
  1733. X
  1734. Xtypedef struct t_spheredata {
  1735. X    Vec    tri_P[3] ;
  1736. X    Vec    tri_N[3] ;
  1737. X    Vec    tri_bb[3] ;
  1738. X} TriData ;
  1739. X
  1740. Xint TriPrint ();
  1741. Xint TriIntersect ();
  1742. Xint TriNormal ();
  1743. X
  1744. XObjectProcs TriProcs = {
  1745. X    TriPrint,
  1746. X    TriIntersect,
  1747. X    TriNormal,
  1748. X} ;
  1749. X
  1750. Xint 
  1751. XTriPrint(obj)
  1752. X Object *obj ;
  1753. X{
  1754. X    TriData *td ;
  1755. X    int i ;
  1756. X
  1757. X    td = (TriData *) obj -> o_data ;
  1758. X    printf("pp 3\n") ;
  1759. X    for (i = 0 ; i < 3 ; i ++) {
  1760. X        VecPrint("\t", td -> tri_P[i]) ;
  1761. X        VecPrint("\t", td -> tri_N[i]) ;
  1762. X    }
  1763. X}
  1764. X
  1765. XTriIntersect(obj, ray, hit)
  1766. X Object * obj ;
  1767. X Ray * ray ;
  1768. X Isect * hit ;
  1769. X{
  1770. X    TriData *td ;
  1771. X    Flt n, d, dist ;
  1772. X    Flt r, s, t ;
  1773. X    Flt a, b ;
  1774. X    Vec P, Q ;
  1775. X
  1776. X    td = (TriData *) obj -> o_data ;
  1777. X
  1778. X    /*
  1779. X     * The matrix td -> tri_bb transforms vectors in the world 
  1780. X     * space into a space with the following properties.
  1781. X     *
  1782. X     * 1.  The sides of the triangle are coincident with the
  1783. X     *     x and y axis, and have unit length.
  1784. X     * 2.  The normal to the triangle is coincident with the 
  1785. X     *     z axis.
  1786. X     *
  1787. X     */
  1788. X
  1789. X    /*
  1790. X     * d is the slope with respect to the z axis.  If d is zero, then
  1791. X     * the ray is parallel to the plane of the polygon, and we count 
  1792. X     * it as a miss...
  1793. X     */
  1794. X
  1795. X    d = VecDot(ray -> D, td -> tri_bb[2]) ;
  1796. X    if (fabs(d) < rayeps)
  1797. X        return 0 ;
  1798. X
  1799. X    /*
  1800. X     * Q is a vector from the eye to the triangles "origin" vertex.
  1801. X     * n is then set to be the distance of the tranformed eyepoint
  1802. X     * to the plane in the polygon.
  1803. X     * Together, n and d allow you to find the distance to the polygon, 
  1804. X     * which is merely n / d.
  1805. X     */
  1806. X
  1807. X    VecSub(td -> tri_P[0], ray -> P, Q) ;
  1808. X
  1809. X    n = VecDot(Q, td -> tri_bb[2]) ;
  1810. X
  1811. X    dist = n / d ;
  1812. X
  1813. X    if (dist < rayeps) 
  1814. X        return 0 ;
  1815. X    
  1816. X    /* 
  1817. X     * Q is the point we hit.  Find its position relative to the
  1818. X     * origin of the triangle.
  1819. X     */
  1820. X
  1821. X    RayPoint(ray, dist, Q) ;
  1822. X    VecSub(Q, td -> tri_P[0], Q) ;
  1823. X
  1824. X    a = VecDot(Q, td -> tri_bb[0]) ;
  1825. X    b = VecDot(Q, td -> tri_bb[1]) ;
  1826. X
  1827. X    if (a < 0.0 || b < 0.0 || a + b > 1.0) 
  1828. X        return 0 ;
  1829. X    
  1830. X    r = 1.0 - a - b ;
  1831. X    s = a ;
  1832. X    t = b ;
  1833. X
  1834. X    hit -> isect_t = dist ;
  1835. X    hit -> isect_enter = 0 ;
  1836. X    hit -> isect_prim = obj ;
  1837. X    hit -> isect_surf = obj -> o_surf ;
  1838. X
  1839. X    VecZero(hit->isect_normal) ;
  1840. X    VecAddS(r, td -> tri_N[0], hit -> isect_normal, hit -> isect_normal) ;
  1841. X    VecAddS(s, td -> tri_N[1], hit -> isect_normal, hit -> isect_normal) ;
  1842. X    VecAddS(t, td -> tri_N[2], hit -> isect_normal, hit -> isect_normal) ;
  1843. X    VecNormalize(hit -> isect_normal) ;
  1844. X
  1845. X    return(1) ;
  1846. X}
  1847. X
  1848. Xint
  1849. XTriNormal(obj, hit, P, N)
  1850. X Object * obj ;
  1851. X Isect * hit ;
  1852. X Point P, N ;
  1853. X{
  1854. X    VecCopy(hit -> isect_normal, N) ;
  1855. X}
  1856. X
  1857. XObject *
  1858. XMakeTri(point)
  1859. X Vec *point ;
  1860. X{
  1861. X    Object * o ;
  1862. X    TriData * td ;
  1863. X    Vec     N, P, Q;
  1864. X    int i, j ;
  1865. X    Flt dmin, dmax, d ;
  1866. X    Vec B[3] ;
  1867. X
  1868. X    o = (Object *) malloc (sizeof(Object)) ;
  1869. X    o -> o_type = T_TRI ;
  1870. X    o -> o_procs = & TriProcs ;
  1871. X    o -> o_surf = CurrentSurface ;
  1872. X
  1873. X    td = (TriData *) malloc (sizeof(TriData)) ;
  1874. X
  1875. X    /* 
  1876. X     * copy in the points....
  1877. X     */
  1878. X    VecCopy(point[0], td -> tri_P[0]) ;
  1879. X    VecCopy(point[2], td -> tri_P[1]) ;
  1880. X    VecCopy(point[4], td -> tri_P[2]) ;
  1881. X
  1882. X    /*
  1883. X     * and the normals, then normalize them...
  1884. X     */
  1885. X    VecCopy(point[1], td -> tri_N[0]) ;
  1886. X    VecCopy(point[3], td -> tri_N[1]) ;
  1887. X    VecCopy(point[5], td -> tri_N[2]) ;
  1888. X    VecNormalize(td -> tri_N[0]) ;
  1889. X    VecNormalize(td -> tri_N[1]) ;
  1890. X    VecNormalize(td -> tri_N[2]) ;
  1891. X
  1892. X    /*
  1893. X     * construct the inverse of the matrix...
  1894. X     * | P1 |
  1895. X     * | P2 |
  1896. X     * | N  |
  1897. X     * and store it in td -> tri_bb[]
  1898. X     */
  1899. X    
  1900. X    VecSub(td -> tri_P[1], td -> tri_P[0], B[0]) ;
  1901. X    VecSub(td -> tri_P[2], td -> tri_P[0], B[1]) ;
  1902. X    VecCross(B[0], B[1], B[2]) ;
  1903. X    VecNormalize(B[2]) ;
  1904. X
  1905. X    InvertMatrix(B, td -> tri_bb) ;
  1906. X
  1907. X    for (i = 0 ; i < NSLABS ; i ++) {
  1908. X        dmin = HUGE ;
  1909. X        dmax = - HUGE ;
  1910. X
  1911. X        for (j = 0 ; j < 3 ; j ++) {
  1912. X            d = VecDot(Slab[i], td -> tri_P[j]) ;
  1913. X            if (d < dmin) dmin = d ;
  1914. X            if (d > dmax) dmax = d ;
  1915. X        }
  1916. X        o -> o_dmin[i] = dmin - rayeps ;
  1917. X        o -> o_dmax[i] = dmax + rayeps ;
  1918. X    }
  1919. X
  1920. X    o -> o_data = (void *) td ;
  1921. X
  1922. X    return o ;
  1923. X}
  1924. X
  1925. XInvertMatrix(in, out)
  1926. X Vec in[3] ;
  1927. X Vec out[3] ;
  1928. X{
  1929. X    int i, j, k ;
  1930. X    Flt tmp, det, sum ;
  1931. X
  1932. X    out[0][0] = (in[1][1] * in[2][2] - in[1][2] * in[2][1]) ;
  1933. X    out[1][0] = -(in[0][1] * in[2][2] - in[0][2] * in[2][1]) ;
  1934. X    out[2][0] = (in[0][1] * in[1][2] - in[0][2] * in[1][1]) ;
  1935. X
  1936. X    out[0][1] = -(in[1][0] * in[2][2] - in[1][2] * in[2][0]) ;
  1937. X    out[1][1] = (in[0][0] * in[2][2] - in[0][2] * in[2][0]) ;
  1938. X    out[2][1] = -(in[0][0] * in[1][2] - in[0][2] * in[1][0]) ;
  1939. X
  1940. X    out[0][2] = (in[1][0] * in[2][1] - in[1][1] * in[2][0]) ;
  1941. X    out[1][2] = -(in[0][0] * in[2][1] - in[0][1] * in[2][0]) ;
  1942. X    out[2][2] = (in[0][0] * in[1][1] - in[0][1] * in[1][0]) ;
  1943. X    
  1944. X    det = 
  1945. X    in[0][0] * in[1][1] * in[2][2] +
  1946. X    in[0][1] * in[1][2] * in[2][0] +
  1947. X    in[0][2] * in[1][0] * in[2][1] -
  1948. X    in[0][2] * in[1][1] * in[2][0] -
  1949. X    in[0][0] * in[1][2] * in[2][1] -
  1950. X    in[0][1] * in[1][0] * in[2][2] ;
  1951. X
  1952. X    det = 1 / det ;
  1953. X
  1954. X    for (i = 0 ; i < 3 ; i ++) {
  1955. X        for (j = 0 ; j < 3 ; j++) {
  1956. X            out[i][j] *= det ;
  1957. X        }
  1958. X    }
  1959. X    
  1960. X#ifdef DEBUG
  1961. X    for (i = 0 ; i < 3 ; i++) {
  1962. X        for (j = 0 ; j < 3 ; j++) {
  1963. X            sum = 0.0 ;
  1964. X            for (k = 0 ; k < 3 ; k++) {
  1965. X                sum += in[i][k] * out[k][j] ;
  1966. X            }
  1967. X            if (fabs(sum) < rayeps) {
  1968. X                sum = 0.0 ;
  1969. X            }
  1970. X            printf(" %g") ;
  1971. X        } 
  1972. X        printf("\n") ;
  1973. X    }
  1974. X    printf("\n") ;
  1975. X#endif /* DEBUG */
  1976. X}
  1977. END_OF_FILE
  1978. if test 6022 -ne `wc -c <'tri.c'`; then
  1979.     echo shar: \"'tri.c'\" unpacked with wrong size!
  1980. fi
  1981. # end of 'tri.c'
  1982. fi
  1983. echo shar: End of archive 2 \(of 3\).
  1984. cp /dev/null ark2isdone
  1985. MISSING=""
  1986. for I in 1 2 3 ; do
  1987.     if test ! -f ark${I}isdone ; then
  1988.     MISSING="${MISSING} ${I}"
  1989.     fi
  1990. done
  1991. if test "${MISSING}" = "" ; then
  1992.     echo You have unpacked all 3 archives.
  1993.     rm -f ark[1-9]isdone
  1994. else
  1995.     echo You still need to unpack the following archives:
  1996.     echo "        " ${MISSING}
  1997. fi
  1998. ##  End of shell archive.
  1999. exit 0
  2000.  
  2001.